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ABSTRACT 


I 
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Several  recent  models  for  the  shock  initiation  of  heterogeneous  explosives 
are  presented,  concentrating  on  those  models  which  have  proved  to  be  the 
most  successful.  Particular  attention  is  given  to  models  of  specific 
interest  to  MRL,  which  are  capable  of  simulating  the  effect  of  particle  size 
on  sensitivity,  and  can  be  readily  incorporated  into  single  phase 
hydrodynamic  computer  codes.  Other  models  are  also  briefly 
considered.  Recommendations  are  made  regarding  the  suitability  of  some 
of  these  models  for  MRL  use. 
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A  CRITICAL  ASSESSMENT  OF  BURN  MODELS  AVAILABLE  FOR 
IMPLEMENTATION  INTO  A  COMPUTER  CODE  TO  MODEL  SHOCK  INITIATION  OF 
HETEROGENEOUS  EXPLOSIVES 


1.  INTRODUCTION 


There  has  been  considerable  interest  for  some  years  at  Materials  Research  Laboratory  in 
the  numerical  simulation  of  shock  and  detonation  phenomena  in  condensed  energetic 
materials.  Our  early  experience  in  non-reactive  shocks  involved  use  of  the  HELP  code,  a 
multimaterial  two-dimensional  Eulerian  code  for  the  numerical  simulation  of  hydrodynamic 
and  elastic-plastic  flow.  HELP  was  used  primarily  for  warhead  modelling  and  was 
particularly  useful  in  the  numerical  simulation  of  the  MRL  38  mm  shaped  charge  jet  [II. 
Whilst  the  general  agreement  between  the  code  calculations  and  the  flash  radiographs  was 
excellent,  problems  were  encountered  in  calculating  an  accurate  shape  for  the  jet  tip 
because  of  the  numerical  scheme  used  for  tracking  the  interface  between  the  different 
materials.  For  these  and  other  reasons,  we  have  now  replaced  HELP  with  the  HULL 
code.  HULL  was  originally  developed  by  Systems,  Science  and  Software  for  the  US  Air 
Force  and  is  now  maintained  by  Orlando  Technology.  The  code  has  both  Eulerian  and 
Lagrangian  capabilities  and  can  perform  simulations  in  either  two  or  three  dimensions. 
HULL  has  an  extensive  library  of  equation-of -state  subroutines  for  non-reactive  media,  but 
only  a  rudimentary  capability  for  numerical  simulation  of  detonation  phenomena. 

Modelling  the  detonation  of  a  condensed  explosive  at  MRL  has  until  now  been 
carried  out  using  the  reactive  hydrocodes  SIN  and  2DL  developed  at  Los  Alamos  National 
Laboratory  by  Mader  in  the  1960s  and  1970s  [2,31.  Both  are  one  and  two-dimensional 
hydrocodes  which  were  specifically  designed  to  model  reactive  shock  phenomena  and  they 
contain  a  variety  of  explosive  bum  subroutines.  However,  significant  advances  have  been 
made  over  the  past  decade  in  the  understanding  of  the  shock  initiation  of  condensed 
heterogeneous  explosives;  SIN  and  2DL  do  not  reflect  these  advances.  In  particular,  the 
most  advanced  bum  subroutine  in  2DL  is  the  Forest  Fire  model  141,  which  is  an  expression 
for  the  rate  of  explosive  decomposition  as  a  function  of  the  local  pressure  which  reproduces 
the  observed  run-to-detonation  distance  for  a  given  initial  input  pressure  (the  so  called  Pop 
plot  data)  (41).  Forest  Fire  has  proved  to  be  very  effective  for  the  numerical  simulation  of 
two  and  three-dimensional  reactive  shock  interactions,  but  cannot  probe  the  basic  physics 
and  chemistry  of  the  shock  initiation  process.  Forest  has  described  a  method  for 
estimating  changes  in  the  Forest  Fire  coefficients  due  to  changes  in  the  explosive  density 
[51,  but  if  the  explosive  particle  size  is  altered  (a  parameter  which  has  an  important 
influence  on  the  shock  initiation  behaviour  of  heterogeneous  explosives  161),  then  a  new 
series  of  experiments  must  be  run  to  generate  new  Pop  plot  data  for  the  Forest  Fire 
coefficients.  Such  experiments  are  very  time  consuming  and  costly,  and  cannot  readily  be 
conducted  at  MRL. 

In  the  last  few  years  a  number  of  models  have  been  proposed  for  the  shock 
initiation  of  heterogeneous  explosives  which  provide  some  insight  into  the  basic  physics 


involved.  One  of  the  most  successful  of  these  is  the  Ignition  and  Growth  model  of  Lee  and 
Tarver  17,81.  In  this  approach  the  initiation  of  detonation  is  divided  into  two  stages:  in  the 
first  stage  the  passage  of  the  shock  front  creates  "hot  spots",  or  localized  areas  of  high 
temperature,  at  density  discontinuities  in  the  heterogeneous  explosive.  In  the  second  stage 
these  hot  spots  grow,  by  a  grain  burning  process,  and  then  eventually  coalesce  to  form  a 
stable  detonation.  While  the  Ignition  and  Growth  model  has  greatly  increased  our 
understanding  of  the  shock  initiation  process,  it  is  still  a  phenomenological  model  in  the 
sense  that  the  coefficients,  which  have  to  be  fitted  to  experimental  embedded  gauge  data, 
need  to  be  redetermined  each  time  one  of  the  physical  properties  of  the  explosive  is 
changed. 


In  a  recent  paper  Kim  has  developed  a  model  which  incorporates  physical 
properties  of  the  explosive  such  as  grain  size  and  porosity  in  a  simple  and  yet  effective 
manner  191.  By  focussing  attention  on  the  basic  physics  of  the  shock  initiation  process  he 
has  been  able  to  develop  a  version  of  the  Ignition  and  Growth  model  which  contains  fewer 
adjustable  constants  than  the  original,  and  which  has  the  advantage  that  once  these 
constants  are  determined  (again,  by  comparison  with  pressure-time  histories  provided  by 
embedded  gauge  data),  they  do  not  have  to  be  changed  to  accommodate  variations  in  grain 
size  or  density/porosity. 

The  purpose  of  this  report  is  to  examine  several  recent  shock  initiation  models 
in  detail,  and  to  assess  their  suitability  for  implementation  in  existing  hydrocodes  at 
MRL.  We  begin  by  first  describing  several  of  the  mechanisms  which  have  been  proposed 
for  the  formation  of  hot  spots  in  heterogeneous  explosives.  Next,  a  critical  examination  of 
the  Ignition  and  Growth  model  of  Lee  and  Tarver  is  made,  followed  by  the  model  of  Kim. 

A  different  approach  is  discussed  in  section  4,  based  on  a  coupled  set  of  rate  equations  for 
the  various  processes  involved,  and  appropriate  time  constants  for  each  of  these 
processes.  In  section  5  other  models  which  have  been  proposed  recently  are  briefly 
described.  The  intention  here  is  not  to  review  all  models  exhaustively,  but  rather  to 
concentrate  on  those  models  of  most  relevance  to  our  own  needs.  In  the  conclusion,  the 
relative  merits  and  disadvantages  of  these  models  are  discussed  and  their  ease  of 
implementation  in  various  hydrocodes  are  considered. 


2.  PROPOSED  MECHANISMS  FOR  HOT  SPOT  FORMATION 


The  shock  initiation  of  a  heterogeneous  explosive  differs  fundamentally  from  that  of  a 
homogeneous  one.  The  differences  were  clearly  described  in  the  early  1960s  in  two 
extensive  experimental  papers  by  the  Los  Alamos  group  of  Campbell,  Davis,  Ramsay  and 
Travis  [6,101.  In  a  homogeneous  explosive,  such  as  a  liquid  or  a  single  crystal,  detonation  is 
caused  by  a  thermal  explosion  due  to  the  bulk  heating  of  the  sample  by  the  passage  of  the 
shock. 


In  a  heterogeneous  explosive  the  situation  is  more  complicated.  Most 
condensed  explosives  consist  of  polycrystalline  materials  containing  voids  of  various  shapes 
and  sizes,  defect  structures,  and  often  small  amounts  of  polymeric  binders  and 
plasticisers.  When  a  shock  wave  travels  through  such  material  it  provides  heating  both  by 
bulk  compression  and  by  the  interaction  of  the  shock  with  the  various  density  discontinuities 
and  defect  structures.  The  localised  regions  of  high  temperature  caused  by  these  density 
discontinuities,  the  hot  spots,  may  then  begin  to  react  and,  if  conditions  are  favourable, 
may  It  ad  to  the  formation  of  a  stable  detonation  even  though  the  temperature  rise  caused 
by  the  bulk  heating  may  be  insufficient  by  itself  to  lead  to  detonation  [171. 
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This  situation  is  illustrated  by  the  shock  initiation  of  the  explosive  PETN  Ill). 
For  a  single  crystal  of  PETN  to  detonate  within  1  us  it  has  to  be  shocked  to  11  GPa,  which 
produces  a  bulk  temperature  of  570 'C.  But  if  the  PETN  is  in  powder  form  and  pressed  to 
almost  crystal  density  it  requires  a  shock  of  only  2.5  GPa.  The  bulk  temperature  in  this 
case  is  only  145" C  and  would  be  insufficient  to  cause  detonation  within  1  us. 

Current  understanding  of  the  initiation  of  detonation  in  heterogeneous 
explosives  by  shock  divides  the  process  into  two  distinct  stages. 

(i)  Ignition  of  a  small  fraction  of  the  explosive  at  random  sites  within  the  sample  due 
to  the  creation  of  hot  spots. 

(ii)  Buildup  to  detonation  from  the  energy  released  by  the  grain  burning  processes 
growing  from  the  original  hot  spots.  If  steady  detonation  is  to  be  attained  it  is 
essential  that  any  losses  in  energy  must  be  less  than  the  energy  released  by  the  hot 
spots. 

There  is  clear  experimental  evidence  supportin  :  this  picture  of  the  initiation  process  [121. 

To  numerically  simulate  shock  initiation  we  therefore  need  to  find  appropriate  models  for 
both  the  ignition  and  buildup  stages.  Many  models  of  hot  spot  formation  have  been 
proposed  for  the  ignition  stage  and  these  are  listed  below. 

1.  Stagnat’on  of  microjets  [131. 

2.  The  hydrodynamic  hot  spot  [141. 

3.  Visco-plastic  heating  in  material  near  the  surface  of  a  collapsing  void  1151. 

4.  Shock  collision  around  high  impedance  inclusions  1161. 

5.  Friction  between  crystal  grains  1171. 

6.  Internal  shear  banding  and  dislocation  pile-up  1181. 

7.  Adiabatic  gas  compression  (171. 

The  buildup  stage  is  usually  modelled  by  a  grain  burning  process  as  first  proposed  by  Eyring 
[191.  There  is  ample  experimental  evidence  for  the  importance  of  this  mechanism  in  the 
buildup  stage  120,211,  and  grain  burning  models  will  be  discussed  in  more  detail  in  the 
following  sections.  In  the  remainder  of  this  section  a  brief  description  is  given  of  some  of 
the  proposed  mechanisms  for  hot  spot  formation  listed  above. 

The  stagnation  of  micro  jets  theory  was  first  proposed  by  Seely  [131.  The  basic 
idea  is  that  as  the  propagating  shock  exits  an  explosive  grain  it  causes  material  to  fly  off 
the  grain  surface.  As  the  grains  are  randomly  oriented  with  respect  to  one  another  the 
material  can  interact  in  various  ways,  and  in  some  cases  jets  will  be  produced.  Seely 
assumes  that  these  jets  will  behave  hydrodynamically,  and  after  traversing  a  typical  void 
space  they  will  collide  with  the  surface  directly  ahead,  and  stagnate.  Seely  uses  simple 
arguments  to  show  that  the  hot  spot  temperature  produced  by  this  mechanism  is 
proportional  to  the  square  of  the  particle  velocity  immediately  behind  the  shock  front.  Lee 
and  Tarver  have  been  able  to  represent  this  mechanism  in  their  original  ignition  and  growth 
model  [71,  and  we  will  comment  on  the  degree  to  which  it  agrees  with  experiment  in  the 
next  section. 

Mader  has  proposed  a  hydrodynamic  void  closure  mechanism  to  explain  hot  spot 
formation  [141.  In  this  model  the  passage  of  the  shock  front  causes  collapse  of  the  void, 
and  high  temperatures  are  produced  by  the  high  impact  pressures  and  focussing  effects 
during  the  collapse  process.  Mader  has  extended  this  model  to  three  dimensions  and  shown 
that  it  is  capable  of  reproducing  the  observed  desensitization  of  heterogeneous  explosives 
by  a  weak  preshock  [221.  In  this  case  the  preshock  closes  the  voids  but  the  hot  spots 
formed  do  not  have  high  enough  temperatures  to  build  up  to  detonation.  The  following 
shock  then  sees  an  effectively  homogenous  material,  and  is  therefore  less  sensitive  to  shock 
initiation. 


9 


The  visco-plastic  heating  of  the  material  in  a  thin  shell  surrounding  the  surface 
of  a  collapsing  void  has  been  considered  by  several  authors  115,23,24).  Frey  1151  has 
presented  quite  a  detailed  description  of  cavity  collapse  in  energetic  materials  based  on 
earlier  work  of  Carroll  and  Holt  1251.  He  has  considered  hot  spot  formation  due  to  inviscid 
plastic  work,  viscoplastic  work,  gas  phase  heating,  and  solid  phase  compression  (ie.  the 
hydrodynamic  model  proposed  by  Mader),  and  presented  a  very  comprehensive  analysis  of 
the  conditions  under  which  each  mechanism  will  occur.  His  main  conclusion  is  that 
viscoplastic  work  is  by  far  the  most  efficient  mechanism  for  producing  high  temperatures, 
and  is  favoured  by  high  viscosity,  low  yield  strength,  and  short  rise  times. 

Frey  was  unable  to  include  the  formation  of  shear  bands  in  his  model  of  cavity 
collapse,  but  has  considered  this  mechanism  in  detail  elsewhere  1261.  Coffey  has  also 
studied  the  effect  of  shear  banding  on  shock  initiation  (27,  281.  These  authors  have  shown 
that  significant  heating  can  occur  in  crystalline  explosives  by  the  generation  and  movement 
of  dislocations  and  the  localization  of  these  dislocations  in  shear  bands  in  the  deforming 
solid.  Direct  experimental  evidence  that  shear  bands  are  important  in  the  initiation  of 
explosives  can  be  found  in  the  paper  by  Field,  Swallowe  and  Heavens  1291.  The  conditions 
under  which  shear  banding  is  the  dominant  mechanism  for  hot  spot  formation  in  any  given 
situation  do  not  seem  to  be  as  well  characterised  as  those  mechanisms  previously  discussed, 
although  it  is  certainly  an  important  mechanism  under  impact  conditions. 

Adiabatic  gas  compression  heating  was  also  considered  by  Frey  1151,  and  was 
shown  to  make  a  negligible  contribution  to  the  temperature  increase  under  typical  shock 
initiation  conditions.  For  impact  conditions  however  adiabatic  gas  compression  can  be 
important. 


3.  THE  IGNITION  AND  GROWTH  MODEL 


In  their  original  model  Lee  and  Tarver  171  divided  the  initiation  process  into  two  distinct 
stages.  In  the  first,  the  ignition  stage,  the  passage  of  the  shock  front  creates  localised 
regions  of  high  temperature  (the  hot  spots)  at  density  discontinuities  within  the 
heterogeneous  material.  In  the  second,  the  growth  (or  build  up)  stage,  these  hot  spots  were 
postulated  to  grow  by  a  grain  burning  process  until  they  eventually  coalesced  to  form  a 
stable  detonation.  The  model  is  phenomenological  in  the  sense  that  plausable  assumptions 
are  made  regarding  the  physical  mechanisms  for  each  of  these  stages  and  then  a  generalized 
energy  release  rate  equation  of  the  form 

=  1(1  -  F>V  +  G(1  -  F)x  F*  Pz  (1) 

Ot 

7  =  PS/PQ  -  1  (2) 


is  considered,  where  F  is  the  fraction  of  explosive  that  has  reacted,  p.  is  the  initial 
density,  pg  is  the  density  of  shocked  explosive,  P  is  the  pressure,  andl,  G,  x,  y,  r  and  z  are 
constants  which  must  be  determined  experimentally. 

Different  models  of  hot  spot  formation  lead  to  different  values  for  the  constant 
r.  In  the  first  reported  application  of  the  model  the  ignition  rate  was  assumed  to  be 
proportional  to  the  strain  rate  in  the  shocked  explosive,  which  corresponds  to  a  value  of  r  = 
1  in  equation  (1).  Some  success  was  achieved  in  modelling  shock  initiation  experiments  in 
PBX-9404  (HMX/NC/TCP  94:3:3),  but  this  value  was  later  found  to  be  innappropriate  over 


10 


* 


* 


an  extended  range  of  stimuli  or  material  properties  171.  A  value  of  r  =  3  corresponds  to  a 
hot-spot  formation  model  due  to  the  stagnation  of  micro  jets,  first  described  by  Seely,  1131, 
while  r  =  4  describes  a  model  based  on  the  amount  of  plastic  work  in  the  void  or  cavity  walls 
required  for  void  collapse.  Best  overall  agreement  with  experiment  was  found  when  r  =  4 
was  used  in  equation  (1)  (71. 

The  constants  x  and  y  are  related  to  the  choice  of  the  geometry  for  the  hot  spot 
combustion  process.  Hot  spots  can  be  considered  to  burn  outward  from  the  void  centre,  or 
inward  over  the  total  grain  surface.  Lee  and  Tarver  considered  a  spherical  hot  spot  burning 
outward,  which  corresponds  to  y  =  2/3.  Requiring  that  the  rate  be  a  maximum  when  the 
combustion  surfaces  overlap  leads  to  the  value  x  =  2/9. 

The  remaining  constants  I,  G  and  z  were  found  by  fitting  to  the  detailed  shapes 
of  embedded  pressure  gauge  records  and  run  distance  to  detonation  data.  Values  of  these 
constants  for  the  three  pressed  explosives  PBX-9404,  TATB,  PETN  and  cast  TNT  are  shown 
in  Table  1.  The  Pz  term  in  equation  (1)  represents  a  laminar  bum  rate  and  the  constant  G 
corresponds  to  a  surface  area  to  volume  ratio.  Measured  values  of  z  for  laminar 
deflagration  rates  in  explosives  are  usually  of  the  order  of  0.8  to  1.0  for  pressures  below 
0.1  GPa.  The  higher  values  for  z  reported  in  Table  1  are  probably  caused  by  grain  fracture 
occuring  due  to  the  higher  pressures  generated.  Lee  and  Tarver  also  used  equation  (1)  with 
a  fixed  value  for  z  of  1.0.  They  found  that  they  could  reproduce  all  the  buildup  and  run 
distance  to  detonation  data  from  embedded  pressure  gauges  by  allowing  the  growth 
coefficient  G  to  increase  as  the  input  shock  pressure  increased. 

Equation  (1)  was  combined  with  the  Jones-Wilkins-Lee  equation  of  state  both  for 
unreacted  explosives  and  their  reaction  products  and  implemented  in  a  one-dimensional 
Lagrangian  hydrodynamic  code.  Using  the  values  of  the  constants  described  above,  and 
shown  in  Table  1,  the  model  was  able  to  successfully  calculate  all  the  sustained  pulse 
manganin  pressure  gauge  and  particle  velocity  data  for  the  four  explosives  at  densities  near 
their  theoretical  maximum  densities.  Lee  and  Tarver  also  found  that  these  same  constants 
gave  a  reasonable  prediction  of  the  shock  initiation  properties  of  the  same  explosives  at 
lower  densities. 

The  model  was  also  applied  to  short  pulse  duration  shock  initiation  experiments 
[71  and  there  it  was  less  successful.  Detailed  quantitative  modelling  of  experimental  data 
was  only  possible  if  the  coefficient  for  the  growth  of  reaction  was  increased  by  a  factor  of 
two  or  three.  To  overcome  these  difficulties  Tarver,  Hallquist  and  Erickson  proposed  a 
three  term  ignition  and  growth  model  181.  Their  expression  for  the  chemical  reaction  rate 
equation  is 

■ff  =  1(1  -  F)b  C2-  -  1  -  a)x  +  G  ( 1  -  F)C  Fd  P* 

at  pQ  1 

+  G2<1  -  F)e  F®  pz .  <31 

The  parameter  a  is  used  to  prohibit  ignition  until  a  certain  degree  of  compression  has  been 
reached.  Values  for  the  constants  I,  Gj,  G2,  a,  b,  c,  d,  e,  g,  x,  y  and  z  for  the  explosives 
PBX-9404  and  LX-17  (TATB/KelF800  92.5:7.5)  are  given  in  Table  2. 

The  main  idea  behind  the  three  term  reaction  rate  model  is  to  split  the  growth 
term  into  two  parts.  The  first  part  models  the  relatively  slow  pressure  and  particle 
velocity  increases  behind  the  shock  front  in  the  early  stages  of  hot-spot  growth,  and 
consequently  the  constant  y  has  the  value  1.0  to  correspond  to  measured  values  in 
deflagration  experiments.  The  second  part  of  the  growth  term  describes  the  rapid 
completion  of  the  reaction  as  the  hot-spots  begin  to  coalesce.  With  Gj  and  G2  fixed,  the 
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value  of  z  is  therefore  greater  than  1.0.  Equation  (3)  was  used  with  the  parameters  given 
in  Table  2  and  yielded  good  overall  agreement  with  a  large  collection  of  shock  initiation 
data  on  PBX-9404  from  three  laboratories  [81. 

While  the  ignition  and  growth  model  has  successfully  calculated  a  great  deal  of 
one  and  two-dimensional  experimenta’  data  on  heterogeneous  explosives  it  is  still  a 
phenomenological  model  in  the  sense  that  it  does  not  attempt  to  model  in  detail  the 
processes  that  cause  the  heterogeneous  heating  of  the  explosive.  Consequently  it  is  not 
capable  of  modelling  the  effect  of  particle  size  or  initial  temperature  on  shock  initiation. 
Cochran  and  Tarver  overcame  this  problem  by  developing  a  statistical  treatment  of  hot¬ 
spot  formation  which  does  take  these  physical  parameters  into  account  and  which  can  be 
used  with  the  ignition  and  growth  model  [301.  They  have  applied  this  combined  model  to 
LX-17  and  have  been  able  to  s\  ccessfully  reproduce  experimental  manganin  pressure  gauge 
records  for  different  particle  sizes,  as  well  as  Pop  plots  determined  at  different  initial 
temperatures. 


4.  THE  KIM  MODEL 


Kim  has  recently  described  an  ignition  and  growth  model  which  incorporates  particle  size 
effects  in  a  natural  manner  [91.  The  model  is  based  on  previous  work  by  Carroll  and  Holt 
[251,  who  showed  that  a  simple  hollow  sphere  model  was  effective  in  describing  the 
compaction  process  in  powdered  materials.  In  the  Kim  model  the  hot-spots  are  assumed  to 
be  generated  around  pores  (voids)  in  the  compacted  explosive.  As  we  have  noted  in 
Section  2,  there  are  many  physical  mechanisms  which  can  lead  to  hot-spot  generation  in 
shocked  porous  explosives.  Kim  and  Soim  [231  have  pointed  out  that,  apart  from  adiabatic 
gas  compression  (which  does  not  contribute  significantly  at  the  high  pressures  short  duration 
shocks  typical  of  shock  initiation),  each  of  these  mechanisms  is  simply  one  form  or  another 
of  mechanical  deformation  of  the  explosive  crystals,  and  in  a  real  explosive  it  is  likely  that 
several  of  these  mechanisms  may  be  operating  simultaneously. 

Faced  with  the  impossibility  of  describing  each  of  these  mechanisms  in  detail, 
Kim  and  Sohn  [231  adopt  a  model  which  is  capable  of  predicting  the  global  behaviour  of  the 
deforming  region.  A  simple  linear  elastic-viscoplastic  relationship  is  assumed  for  the 
material  behaviour: 


de 

dt 


y^a 


.  1  do 

V  +  E  dt 


(4) 


where  e  is  the  stra'n,  o  the  stress,  and  o  the  static  yield  strength  of  the  material.  The 
first  term  on  the  right  hand  side  of  equation  (4)  describes  the  viscoplastic  response, 
with  y  being  the  coefficient  of  viscosity.  The  second  term  represents  the  elastic 
response,  and  E  is  the  Young’s  modulus  of  the  material. 

The  model  adopted  for  the  compaction  process  assumes  that  pores  within  the 
explosive  are  uniform  in  size  and  distribution.  Figure  1  illustrates  how  a  typical  void  in  a 
porous  explosive  is  represented  in  the  hollow  sphere  model.  If  the  average  particle 
diameter  d  has  a  value  of  200  jxn,  then  the  radius  of  the  outer  sphere  rQ  (also  shown  as  a 
dotted  line-in  Fig.  la)  is  chosen  to  be  100  ttn.  If  the  explosive  has  a  porosity  of  2%  then 
the  value  of  the  inner  radius  is  determined  by  requiring  that  the  ratio  of  the  volume  of  the 
inner  sphere  to  that  of  the  outer  sphere  should  be  0.02,  leading  to  a  value  for  r-  of  (0.02)1  3 
r_.  The  equations  of  motion  for  the  compression  of  the  material  surrounding  the  pore  have 
the  form 
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(5) 


do 
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+  2 


=  0 


dv 

dr 


i  ■ 
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l  i_  , 

2G  at  (or 


(6) 


dv 
3  r 


<7> 


where  a  is  the  radial  stress,  o0  is  the  tangential  stress,  r  is  the  space  coordinate,  v  is  the 
radial  displacement  rate, 


y  -  k  =  Oq / \/3 ,  and  G  =  E/3. 


These  can  be  solved  analytically  using  Laplace  transforms,  and  an  expression  for  the 
temperature  rise  within  the  material  due  to  bulk  mechanical  deformation  can  be  written  as 
follows: 


dT 

dt 


<P„  -  P  -  2>/3  kun  ~r> 

0  g  r.  V 


,  -3  -3,2  6 

(ri  "  r0  ’  r 


_1 
'  i/3 


(8) 


Here  p  is  the  solid  density,  C  its  heat  capacity,  P_  is  the  applied  stress  at  rQ  and  P  is  the 
gas  pressure  in  the  void.  This  temperature  rise  within  the  explosive  due  to  this  mechanical 
deformation  leads  to  both  heat  conduction  and  chemical  reaction.  The  equation  describing 
the  overall  temperature  change  can  therefore  be  written  as 


o  C 


3T(  r , t ) 

at 


<S> 


a_ 

3r 


(r2  k 


*  JT 
ar 


)  +  Q 


a  a  (r.t) 
at 


(9) 


where  k*  is  the  thermal  conductivity,  Q  is  the  heat  of  reaction,  and  t  is  the  degree  of 
reaction  of  the  explosive,  with  A  =  1  when  the  explosive  is  fully  reacted.  The  last  term  in 
equation  (9)  is  taken  to  have  a  simple  Arrhenius  reaction  rate 


3A( r ,  t  > 

at 


(1  -  A)  Z  exp 


* 


T(r,t) 


(10) 


where  Z  is  the  frequency  factor  and  T*  is  the  activation  temperature. 

The  model  under  discussion  here  is  a  microscopic  model  in  the  sense 
that  A (  r ,  t )  represents  a  local  reaction  rate  within  the  shell  of  material  surrounding  the 
hollow  sphere,  showing  a  maximum  at  the  inner  radius  r-  and  a  minimum  at  the  outer  radius 
rg,  and  a  monotonic  decrease  between.  Use  of  this  expression  in  a  macroscopic 
hydrodynamic  code  requires  it  to  be  first  averaged  over  the  microscopic  coordinate  r.  The 
local  reaction  rate  is  hence  integrated  over  the  material  within  the  shell  to  define  a 
macroscopic  reaction  rate  a«  fu.iows 
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j* 

A(t).  =  f  0  A  (r,t)  47tt2  dr  /  4  (r  3  -  r.^).  (11) 

ign  Jr.  3  o  i 


We  now  consider  the  terms  used  to  describe  the  growth  of  the  reaction  after 
ignition  by  the  mechanical  deformation  described  above.  In  the  early  stages  of  chemical 
reaction  it  is  assumed  that  the  reaction  progresses  outward  from  the  surface  of  the  hollow 
sphere  in  a  surface  burning  mode.  The  reaction  rate  can  therefore  be  expressed  as 


dA  =  AF^  A2/3 
dt  r 

o 


(12) 


where  A  and  y  are  constants  which  describe  the  slow  burning  rate  at  the  low  pressures  found 
at  the  beginning  stage  of  the  reaction.  After  the  reaction  has  progressed  to  a  significant 
degree  the  burning  will  have  penetrated  all  the  crevices  between  the  particles  and  the 
reaction  will  then  be  dominated  by  inward  burning  of  the  individual  particles.  This  change 
in  geometry  results  in  a  reaction  rate  of  the  form 


dA  _  BP2  (1  -  A)2/3 
dt  r 


where  now  the  constants  B  and  z  describe  the  fast  burning  rate  at  the  high  pressures 
generated  near  the  completion  of  the  reaction. 

Only  five  coefficients  need  to  be  fitted  to  experimental  data  to  describe  the 
shock  initiation  of  a  porous  explosive  using  this  model.  These  are  the  constants  A,  B,  y,  z, 
and  the  material  viscosity  parameter  y1  (the  constant  E  was  found  to  have  no  effect  on  the 
results  and  is  no  longer  considered);  air  other  constants  are  known.  Kim  has  applied  the 
model  to  shock  initiation  data  for  PBX-9404  and  found  good  agreement  with  pressure-time 
histories  obtained  from  embedded  stress  gauges.  He  has  also  considered  an  imaginary 
explosive  with  a  particle  size  four  times  smaller  than  in  the  PBX-9404  study  and  explicitly 
shown  that  the  explosive  with  the  smaller  particle  size  is  harder  to  initiate,  but  transits  to 
detonation  more  readily,  than  the  explosive  with  the  larger  particle  size.  This  behaviour  is 
in  agreement  with  known  experimental  results  [461.  Calculations  are  currently  underway  to 
model  particle  size  effects  in  the  explosive  HNS  131). 


5.  THE  EMPIRICAL  HOT  SPOT  MODEL 


The  empirical  hot  spot  model  has  been  developed  in  a  series  of  papers  by  Johnson,  Tang  and 
Forest  of  Los  Alamos  National  Laboratory  [32-361.  Unlike  the  approach  of  Kim,  and  to  a 
lesser  extent  that  of  Lee  and  Tarver,  no  attempt  is  made  to  explicitly  model  any  of  the 
physical  mechanisms  leading  to  hot-spot  formation.  Instead,  time  constants  are  defined  for 
each  of  the  significant  processes  involved  in  the  overall  shock-to-detonation  transition,  and 
a  series  of  rate  equations  written  down  which  govern  the  time  depsndence  of  each  of  these 
processes.  This  is  achieved  using  the  following  approach. 

Let  t  represent  the  characteristic  time  for  hot  spot  excitation  due  to  the 
passage  of  the  shlick  front.  The  decomposition  of  these  hot  spots  is  then  characterized  by 
the  time  constant  t  .  The  process  for  the  transport  of  energy  from  the  hot  products  of 
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the  hot  spot  decomposition  to  the  colder  bulk  explosive  is  then  governed  by  the 
characteristic  time  t  .  Finally,  is  the  characteristic  time  for  decomposition  in  the 
bulk  explosive  once  tfie  transfer  of  energy  from  the  hot  spots  has  taken  place.  A  set  of  six 
coupled  rate  equations  can  then  be  written  down  to  describe  the  time  dependence  of  the 
reactants,  intermediates,  and  products,  for  both  hot  spots  and  the  remainder  of  the 
explosive.  This  constitutes  a  very  general  model  for  the  overall  shock-to-detonation 
transition. 

By  making  some  plausable  physical  assumptions  Tang,  Johnson  and  Forest  1331 
are  able  to  considerably  simplify  the  general  model.  They  assume  that  the  shock  process 
which  leads  to  the  formation  of  the  hot  spots  will  be  very  much  faster  than  the  process 
leading  to  their  decomposition,  ie.  they  make  the  approximation 


Also,  they  assume  that  the  energy  transfer  process  from  the  hot  products  of  the  hot-spot 
decomposition  to  the  balance  of  the  explosive  will  be  much  slower  than  the  decomposition 
of  the  remainder  of  the  explosive,  ie.  they  make  the  approximation 


t  »  r ,  (15) 

m  a 

With  the  use  of  conditions  (14)  and  (15)  the  six  coupled  rate  equations  reduce  to  just  two 
coupled  equations,  these  being 


—  (1  -  A.) 
rc  h 


dA, 

Hhrr>  (17) 

m  0 

where  p  is  the  hot  spot  mass  fraction,  fg  is  a  threshold  value  of  the  normalized  mass 
fraction  of  hot  spot  reaction  that  must  be  reached  before  the  bum  can  propagate  into  the 
explosive,  A,  is  the  mass  fraction  of  hot  spot  products  divided  by  p ,  and  A,  is  the  mass 
fraction  of  tne  products  of  the  balance  of  the  explosive  divided  by  (1  -  p) ,  The  overall 
reaction  rate  variable  A  is  given  by 


A  =  pA^  +  (1  -  p)  A^, 


and  its  rate  of  change  governed  by  the  equation 


=  -^  (1  -  A,  )  +  ^  (d  -  A)  -  p  (1  -  A,  )]  i 


dt  t 


.  h  *0, 


Use  of  equations  (16)  and  (19)  requires  that  expressions  must  be  found  for  the 
characteristic  times  r  and  t  .  In  the  empirical  hot-spot  model  it  is  assumed  that  the 
passage  of  the  initial  Aock  wave  of  pressure  amplitude  P.  produces  an  average  hot-spot 
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temperature  Tg  given  by 


T0  P8  -1 

Ts  =  To  11  -  m  .  in  <P • 
s  u  T  O 


(20) 


where  m,  Tq  and  Pq  are  constants  and  T*  is  the  Arrhenius  activation  temperature.  The 
characteristic  time  c  is  then  identified  with  the  induction  time  for  thermal  explosion  of 
the  hot  spot,  given  by 


z 
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s  ,T  , 

exp  (=-) , 

T  0Z  s 


(21) 


where  $  is  the  temperature  coefficient  due  to  chemical  reaction  and  Z  is  the  frequency 
factor  for  Arrhenius  reaction. 


The  time  z  characterizes  the  energy  transfer  from  the  hot-spot  products  to  the 
bulk  of  the  unbumed  explosive.  At  low  pressures  it  is  assumed  that  this  process  will  occur 
by  normal  heat  conduction,  while  at  higher  pressures  it  is  expected  that  both  convection  and 
turbulence  could  become  important.  The  expression  for  r  is  therefore  composed  of  two 
parts:  m 


r  =  1G_  P  +  G  (P)  1  1 .  (22) 

m  V 

The  linear  term  represents  the  weaker  energy  transfer  mechanism  while  G(P)  represents  the 
mechanism  dominant  at  higher  pressures.  Tang,  Johnson  and  Forest  134]  in  fact  show  that 
G(P)  should  be  identified  with  the  Forest  Fire  rate,  and  so  has  the  form 

n 

G(P)  =  exp  (  l  a.P1 ) .  (23) 

i=0 

The  empirical  hot-spot  model,  defined  by  equations  (16)  and  (19-23),  has  been 
applied  to  the  explosive  PBX-9404  and  calibrated  to  time-resolved  particle-velocity  data 
for  sustained  shocks  132).  A  one-dimensional  reactive  characteristics  code  was  used  for  the 
calibration,  and  the  model  was  able  to  successfully  reproduce  all  available  particle-velocity 
and  Pop  plot  data.  The  same  model  was  then  used  (with  identical  rate  parameters)  to 
simulate  the  effect  of  a  finite-duration  pulse  in  PBX-9404  [331,  and  good  agreement  with 
experiment  was  again  found.  The  model  has  also  been  implemented  in  the  two-dimensional 
hydrodynamic  finite-element  code  DYNA2D  and  used  to  investigate  the  two-dimensional 
effects  of  comer  turning  and  shock  desensitization  in  PBX-9502  (TATB/Kel-F800  95:5)  and 
again  the  results  have  agreed  well  with  experiment  1331. 

In  a  more  recent  application  of  the  model  Tang  et  al  1341  have  shown  that  it  is 
capable  of  reproducing  the  effect  of  changes  in  both  density  and  grain  size.  Arguments 
were  presented  to  show  that  initiation  behaviour  at  lower  density  could  be  simulated  by 
increasing  the  value  of  the  hot-spot  mass  fraction  u  and  the  reference  temperature  T.  . 
Excellent  agreement  with  experimental  Pop  plot  data  for  PBX-9404  at  two  different  0 
densities  was  obtained  using  this  approach.  The  effect  of  a  change  in  grain  size  on  the 
fundamental  parameters  of  the  theory  is  more  subtle.  Tang  et  al  argue  that  decreasing 
grain  size  should  lead  to  an  increase  in  both  the  sensitivity  parameter  m  and  the  hot-spot 
mass  fraction  j j,  but  decrease  in  the  hot-spot  reference  temperature  T.  .  Changing  these 
constants  in  the  directions  indicated  does  lead  to  excellent  agreement  Between  experiment 
and  simulation  for  Pop  plot  data  on  superfine  and  micronized  TATB.  Tang  has  also  used  the 
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model  to  simulate  grain  size  effects  on  the  performance  of  boosters  [351,  and  extended  the 
model  to  include  a  slow  process  occuring  near  the  end  of  the  reaction  [361. 


6.  OTHER  MODELS 


In  this  section  we  discuss  other  approaches  to  the  shock  initiation  of  heterogeneous 
explosives  which  have  been  described  recently.  A  straight  forward  approach  to  the 
determination  of  the  reaction  rate  law  in  energetic  materials  has  been  described  by 
Anderson  et  al  (371.  This  method  is  based  on  a  Lagrangian  analysis  of  experimental  data 
obtained  from  embedded  manganin  pressure  gauges.  The  rate  law  found  using  this 
empirical  approach  has  come  to  be  known  as  the  DAGMAR  model  (for  Direct  Analysis 
Generated  Modified  Arrhenius  Rate).  It  has  the  form 

ft  =  (1  -  A)  z0  psn  exp  (-TV  (24) 

where  A  is  the  fraction  of  material  reacted,  P  the  shock  strength,  T  is  the  temperature, 
and  Zq,  n,  and  T*  are  constants.  Equation  (24)  has  been  found  to  provide  a  good  description 
of  the  initiation  behaviour  of  PBX-9404  [381  and  porous  TATB  1371.  Note  that  the  pressure 
term  appearing  in  equation  (24)  is  the  initial  shock  strength,  rather  than  the  current 
pressure  P.  Anderson  et  al  found  inclusion  of  a  shock  strength  term  was  crucial  for 
accurate  simulation  of  short-shock  experiments.  While  equation  (24)  is  a  useful  description 
of  the  rate  law  for  the  explosives  mentioned,  it  does  not  describe  grain  size  effects,  nor  the 
mechanism  by  which  the  hot  spots  are  initially  created. 

Another  expression  for  the  reaction  rate  law  which  includes  a  dependence  on  the 
shock  strength  has  been  described  by  Damamme  and  Missonier  [391  and  is  known  as  the 
Krakatoa  model.  It  has  the  form 

=  (NQ)1/3  K(P)  H(  A) ,  (25) 

where  Nq  is  the  number  of  hot  spots  per  unit  volume,  K  is  the  radial  speed  of  growth  of  the 
decomposing  explosive  zone,  and  H  is  a  function  which  depends  on  the  shape  of  the  hot 
spots.  The  model  uses  the  following  functional  forms: 

(N„)1/3  =  A  exp  (I/I  )  (26) 

0  m  ^  a 

where  I  is  a  function  which  depends  on  the  shock  strength,  and  Am  and  Ig  are  constants. 


K(P)  =  Pn 


(27) 


and 


9  /*t 

H(  A)  =  (1  -  A)  I  J>n  (1  -  A>r  . 


(28) 
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The  Krakatoa  model  has  been  applied  to  a  TATB  based  composition  and  by  suitable  choice 
of  the  constants  A,  Ia  and  n  has  been  able  to  successfully  reproduce  experimental  Pop  plot 
data. 

All  the  models  for  shock  initiation  of  heterogeneous  explosives  discussed  so  far 
have  been  essentially  single  fluid  models,  i.e.  the  reacting  material  has  been  treated  as  a 
mixture  of  two  coexisting  phases;  a  solid  phase  consisting  of  unbumed  explosive,  and  a 
gaseous  phase  consisting  of  reaction  products.  A  single  reaction  variable  describes  the 
amount  of  each  phase  present,  and  then  the  specific  volume  and  energy  is  defined  as  a 
weighted  sum  of  these  variables  in  the  separate  phases.  A  suitable  equation  of  state  is 
used  for  each  phase  and  then  pressure  and  temperature  equilibrium  between  the  two  phases 
is  assumed.  Kip  et  al  have  adopted  quite  a  different  approach  1401.  They  have  developed  a 
model  based  on  the  continuum  theory  of  chemically  reacting,  multiphase  mixtures,  and 
applied  it  to  the  heterogeneous  explosive  PBX-9404.  Their  model  for  this  particular 
explosive  consists  of  three  phases:  phase  one  is  the  HMX  granules,  phase  two  the 
nitrocellulose  binder,  and  phase  three  the  reaction  products.  The  TCP  is  simply  treated  as 
an  inert  additive.  Each  phase  has  its  own  pressure,  density,  temperature  and  entropy,  and 
the  volume  fractions  ip .  of  each  phase  are  treated  as  independent  variables.  Separate 
equations  for  conservation  of  mass,  momentum  and  energy  are  solved  for  each  phase,  and 
various  functions  C.  describing  mass  and  energy  transfer  between  phases  must  be  defined. 

In  the  application  to  PBX-9404  it  is  assumed  that  the  hot-spots  form  in  the  binder,  which  is 
assumed  to  break  down  according  to  Arrhenius  kinetics,  ie. 

C2+  =  -  Z2  <p2k2  exp  (-  T*/T2),  (29) 

where  k  is  the  binder  density,  Z2  is  the  frequency  factor,  and  T*  the  activation 
temperature.  The  HMX  particles  decompose  according  to  a  grain  burning  model,  and  so  the 
mass  transfer  rate  for  this  phase  is  given  by 

Cx+  =  -  W1(^32/3/r°)  4,1k1  P3n,  (30) 

where  r®  is  the  average  particle  radius,  is  the  bum  velocity,  and  the  pressure  exponent  n 
is  taken  to  be  a  function  of  the  temperature  of  the  HMX  particles.  This  model  was 
implemented  in  a  one-dimensional  Lagrangian  finite  difference  code  and  numerical 
predictions  were  compared  with  a  substantial  collection  of  shock  and  ramp-wave  data. 

Good  overall  agreement  between  the  simulations  and  experimental  results  was  found.  It 
should  be  noted  here  that  the  multiphase  model  for  PBX-9404  just  described  is  somewhat 
atypical  in  that  most  explosives  do  not  have  an  energetic  binder  such  as  nitrocellulose. 

Baer  and  Nunziato  [411  have  also  used  a  two-phase  model  similar  to  the  one  just 
described  to  model  the  deflagration-to-detonation  transition  in  one  dimensional  columns  of 
HMX.  Khasainov  et  al  [241  have  also  used  a  two-phase  model  to  simulate  the  shock 
initiation  behaviour  of  porous  PETN.  Hot-spot  formulation  is  modelled  by  visco-plastic 
heating  around  individual  pores,  and  growth  of  the  reaction  occurs  by  surface  grain 
burning.  The  properties  of  the  solid  and  gaseous  phases  are  quite  different,  and  are 
modelled  using  a  two-phase  flow  theory.  Apart  from  the  two-phase  formulation,  this  is 
similar  in  some  ways  to  the  model  described  by  Kim  [91.  The  recent  work  of  Nutt  [421  is 
also  relevant  here  as  well.  He  has  also  modelled  hot-spot  formation  by  viscous  heating 
around  pore  collapse,  and  has  incorporated  thiB  with  the  original  growth  model  of  Lee  and 
Tarver  [71,  and  a  detailed  thermodynamic  analysis  of  a  two  component  reactive  system. 
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7.  DISCUSSION 


MRL  interest  in  shock  initiation  is  currently  centered  on  surface  area  effects  in  the  porous 
explosive  HNS  for  use  in  slapper  detonators,  and  particle  size  effects  in  RDX  for  use  in 
PBXs  and  booster  explosives.  Whichever  model  is  adopted  to  describe  these  effects  must 
be  capable  of  implementation  in  an  existing  hydrocode  at  MRL.  As  previously  discussed, 
only  the  reactive  hydrocodes  SIN  and  2DL  have  been  used  to  model  the  behaviour  of 
explosives.  Both  are  single  phase  codes  which  model  both  the  solid  phase  reactants  and  gas 
phase  products  using  the  HOM  equation  of  state  [41.  More  advanced  hydrocodes  are  being 
developed  for  reactive  flow  based  on  the  Flux  Corrected  Transport  technique  144,471,  but 
initial  developments  will  result  in  single  phase  codes.  Because  of  this  limitation,  shock 
initiation  models  which  can  be  implemented  in  single  phase  codes  are  of  particular  interest. 

The  critical  assessment  detailed  in  previous  sections  has  shown  that  there 
currently  exist  a  variety  of  different  models  of  varying  degrees  of  sophistication  for  the 
numerical  simulation  of  shock  initiation  of  energetic  materials.  Because  it  is  virtually 
impossible,  with  the  current  level  of  understanding,  to  clearly  identify  which  of  the  many 
possible  mechanisms  will  be  operable  in  the  ignition  and  growth  processes  in  a  given 
explosive  and  or  mode  of  initiation,  each  model  requires  a  number  of  adjustable 
constants.  The  values  for  each  of  these  constants  are  found  by  fitting  to  experimental 
dynamic  data,  usually  pressure-time  plots  or  particle  velocity-time  plots  at  specific  stations 
within  the  explosive.  This  is  an  iterative,  time  consuming  process,  and  so  most  models 
have  been  fitted  to  only  a  few  explosives. 

The  Ignition  and  Growth  model  [7,81  has  probably  been  fitted  to  more  explosives 
than  any  other  model.  However,  it  probably  also  requires  the  determination  of  more 
adjustable  constants  than  any  other  model.  For  example,  the  three  term  version  described 
by  Tarver,  Hallquist  and  Erickson  requires  the  specification  of  twelve  constants  [81,  and 
more  would  be  required  if  particle  size  effects  were  to  be  modelled  using  the  approach  of 
Cochran  and  Tarver  (301. 

The  Kim  model  is  similar  to  the  Ignition  and  Growth  model  but  contains  a  more 
effective  formulation  of  both  the  ignition  and  growth  terms,  and  consequently  requires 
specification  of  only  five  adjustable  constants.  It  is,  however,  more  computationally 
intensive  than  the  Ignition  and  Growth  model  as  the  heat  equation  has  to  be  solved  within 
each  hydrodynamic  cell  for  each  time  step. 

The  Empirical  Hot  Spot  model  has  been  extensively  developed  and  successfully 
applied  by  its  originators  to  a  range  of  explosives  and  experimental  situations  [32-361.  It 
has  successfully  modelled  both  grain  size  and  density  effects  (341,  but  its  application  in  this 
area  is  not  straightforward.  Subtle  arguements  are  needed  to  indicate  in  which  direction 
various  parameters  should  be  changed,  and  fine  tuning  of  the  parameter  values  is  required 
to  fit  the  experimental  results. 

The  multiphase  models  have  also  been  extensively  developed  and  applied  to  a 
variety  of  explosives  [24,40,41,42,1.  These  models  are  considerably  more  detailed  than  the 
single  phase  models  just  described,  and  consequently  are  capable  of  providing  a  better 
description  of  the  various  processes  occurring  during  shock  initiation.  They  also  require 
specification  of  many  more  physical  parameters,  and  for  this  reason  are  susceptible  to 
problems  of  non-uniqueness  [451.  The  numerical  coding  for  multiphase  models  is  also 
considerably  more  complicated  than  for  the  single  phase  models. 

Because  of  MRL  interest  in  modelling  surface  area  and  particle  size  effects, 
both  the  Ignition  and  Growth  model  and  the  Hot  Spot  model  are  inappropriate.  The  Ignition 
and  Growth  model  was  not  specifically  designed  to  model  particle  size  effects.  These  can 
be  included,  using  the  method  of  Cochran  and  Tarver  [301,  but  by  then  the  number  of 
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adjustable  constants  becomes  unwieldy,  and  are  too  numerous  to  be  experimentally 
determined  at  MRL.  The  Hot  Spot  model  has  been  used  successfully  to  simulate  particle 
size  effects  (341,  but  the  process  is  not  straightforward,  and  requires  a  high  degree  of 
familiarity  with  the  model.  Both  the  DAGMAR  and  Krakatoa  models  are  also 
inappropriate,  as  neither  addresses  the  problem  of  particle  size  effects.  The  multiphase 
models  briefly  discussed  in  Section  6  must  also  be  ruled  out  because  we  currently  do  not 
have  the  multiphase  code  capability  to  implement  them. 

The  most  promising  approach  is  the  Kim  model.  This  has  the  advantage  of 
having  been  specifically  formulated  to  model  particle  size  effects,  as  well  as  requiring  the 
experimental  determination  of  only  a  small  number  of  adjustable  constants.  It  also  uses 
the  HOM  equation  of  state  for  both  reactants  and  products,  which  is  the  equation  of  state 
used  by  both  the  SIN  and  2DL  codes.  The  disadvantage  with  this  model,  as  previously 
noted,  is  that  the  heat  conduction  equation  has  to  be  solved  within  the  shell  of  material 
surrounding  the  hollow  sphere  within  each  hydrodynamic  cell  for  each  hydrodynamic  time 
step.  This  is  not  a  trivial  consideration,  and  makes  the  implementation  of  the  model  into 
an  existing  hydrocode  far  from  straightforward  (311. 

At  MRL,  the  Kim  model  could  currently  be  implemented  into  the  ID  Lagrangian 
finite  difference  code  SIN.  This  code  has  the  advantage  of  being  well  documented  [41,  but 
has  the  disadvantage  that  it  uses  the  artificial  viscosity  method  to  handle  shock  formation 
(431.  This  successfully  dampens  oscillations  around  the  shock  front,  but  can  also  lead  to 
problems  when  parameters  reflecting  real  viscosity  effects  need  to  be  fitted  to 
experimental  data.  Models  of  hot  spot  formation  based  on  viscoplastic  heating  are  a  case 
in  point,  and  problems  of  this  type  have  already  been  encountered  with  the  Kim  model  [311 
during  its  current  implementation  in  a  ID  finite  difference  Lagrangian  code  known  as 
STEALTH  [231,  which  also  uses  the  artificial  viscosity  method  for  shock  formation. 

These  problems  could  be  overcome  by  implementing  the  Kim  model  in  a  ID  Flux 
Corrected  Transport  (FCT)  code  which  is  currently  being  designed  at  MRL  [441.  FCT 
algorithms  require  no  artificial  viscosity  terms  and  the  shock  front  can  be  infinitely 
sharp.  This  would  make  identification  of  an  accurate  value  for  the  material  viscosity  in 
the  Kim  model  much  easier. 

Most  hydrocodes  used  by  the  explosives  community  to  date  use  the  artificial 
viscosity  method  to  handle  shock  formation.  However,  there  have  been  many 
improvements  in  recent  years  in  the  numerical  schemes  developed  for  solving  the  equations 
of  ideal  compressible  flows,  and  none  of  these  techniques  requires  the  use  of  an  artificial 
viscosity  term.  The  explosives  community  has,  in  general,  been  slow  to  adapt  these  newer 
methods  to  the  more  complex  reacting  flows.  One  exception  is  the  work  of  Oran,  Boris  and 
their  collaborators  1481.  They  have  used  the  finer  shock  resolution  available  with  FCT 
algorithms  to  gain  valuable  insights  into  the  stability  of  detonation  in  gaseous  and  liquid 
explosives.  Two  papers  in  the  coming  Ninth  Symposium  (International)  on  Detonation 
[49,501  also  show  the  gradual  adoption  of  these  newer  techniques.  One  of  the  advantages  of 
these  newer  methods  is  that  the  finer  shock  resolution  allows  a  more  accurate 
determination  of  the  radius  of  curvature  of  a  two  dimensional  detonation  front.  Recent 
theories  connecting  curvature  with  reaction  rate  then  allow  information  to  be  obtained  on 
the  rate  of  energy  release  within  the  explosive  [50,511.  Information  of  this  type  is 
currently  of  interest  to  MRL  for  some  recent  PBX  formulations. 

In  conclusion,  the  Kim  model  should  be  used  to  simulate  particle  size  effects, 
and  the  model  should  be  implemented  in  either  the  SIN  hydrocode,  or  an  FCT  hydrocode 
currently  under  development  at  MRL.  The  FCT  hydrocode  is  preferred,  for  reasons 
discussed  above,  but  implementation  in  either  hydrocode  will  not  be  straightforward. 
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Table  1  Values  of  constants  for  original  Ignition  and  Growth  model  (71. 


EXPLOSIVE 

PBX-9404 

TATB 

PETN 

CAST  TNT 

I  ( Is) 

44 

50 

20 

50 

G  ( Is  Mbars 

200 

125 

400 

40 

Z 

1.6 

2.0 

1.4 

1.2 

Table  2  Values  of  constants  for  modified  Ignition  and  Growth  model  18). 


EXPLOSIVE 

PBX-9404 

LX-17 

I  (IS) 

7.43  x  1011 

4.0  x  106 

b 

0.667 

0.667 

a 

0.0 

0.22 

X 

20.0 

7.0 

Cj  (Mbars  y  is ) 

3.1 

0.6 

c 

0.667 

0.667 

d 

0.111 

0.111 

y 

1.0 

1.0 

G2  (Mbars  z  is ) 

400.0 

400.0 

e 

0.333 

0.333 

e 

1.0 

1.0 

Z 

2.0 

3.0 

TABLE  OF  SYMBOLS 


F 

P 


n 

I,  G,  r, 
x,  y,  z 


G-i,  Gn, 
a,  b,  c,  d, 
e,  g 


£ 


O 


Fraction  of  reacted  explosive  (Ignition  and  Growth  model) 

Density,  p  is  initial  density,  p  is  density  of  shocked 
explosive. 


Constants  in  the  Ignition  and  Growth  model  which  must  be 
determined  by  fitting  to  experiment. 

Constants  in  the  modified  Ignition  and  Growth  model  which 
must  be  determined  by  fitting  to  experiment. 


Strain 

Stress,  a  is  the  static  yield  strength,  a  radial 
stress,  o  tangential  stress. 

Coefficient  of  viscosity 

Young’s  modulus 

Radial  space  coordinate 

Radial  velocity 

yl°o 

°oNZ 

Thermal  conductivity 

Inner  and  outer  radii  of  hollow  sphere  used  in  Kibong  Kim 
model. 

Pressure.  PQ  is  applied  stress  at  rQ,  P„  is  gas  pressure 
in  the  void.  6 

Heat  capacity  at  constant  pressure. 

Heat  of  reaction 

Temperature.  T*  is  an  activation  temperature. 

Reaction  rate  variable 

Frequency  factor  in  Arrhenius  reaction  rate. 
Characteristic  time  for  hot  spot  excitation. 
Characteristic  time  for  hot  spot  decomposition 


TABLE  OF  SYMBOLS 

(continued) 


Characteristic  time  for  transport  of  energy  from  hot  spots 
to  colder  bulk  explosive. 

Characteristic  time  for  decomposition  of  bulk  explosive. 

Hot  spot  mass  fraction 

Threshold  value  of  the  normalized  hot  spot  mass  fraction. 

Mass  fraction  of  hot  spot  products  divided  by  u  ■ 

Mass  fraction  of  the  products  of  the  balance  of  the  explosive 
divided  by  (X  -  . 

Shock  strength 

Constant  in  Arrhenius  expression 
Number  of  hot  spots  per  unit  volume 
Volume  fraction  of  phase  i 
Mass  transfer  rate  for  phase  i 
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(b)  Idealization  in  Hollow  Sphere  Model. 


Figure  1  The  Hollow  Sphere  Model 
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